Visualizing raw trial looking time

d_no_target <- d %>% 
  filter(trial_stimulus_type != "target")

Looking time distribution

RAW

data_with_demog %>% 
  ggplot(aes(x = trial_looking_time)) + 
  geom_histogram(bins = 90) + 
  xlim(0, 9000)
## Warning: Removed 72 rows containing non-finite values (stat_bin).
## Warning: Removed 2 rows containing missing values (geom_bar).

data_with_demog %>% 
  ggplot(aes(x = trial_looking_time)) + 
  geom_histogram(bins = 90) + 
  xlim(0, 9000) + 
  facet_wrap(~subject)
## Warning: Removed 72 rows containing non-finite values (stat_bin).
## Warning: Removed 134 rows containing missing values (geom_bar).

AFTER EXCLUSION

d %>% 
  ggplot(aes(x = trial_looking_time)) + 
  geom_histogram(bins = 90) + 
  xlim(0, 6000)
## Warning: Removed 2 rows containing missing values (geom_bar).

d %>% 
  ggplot(aes(x = trial_looking_time)) + 
  geom_histogram(bins = 90) + 
  xlim(0, 6000) + 
  facet_wrap(~subject)
## Warning: Removed 84 rows containing missing values (geom_bar).

Exclude the target trial

d_no_target %>% 
  ggplot(aes(x = trial_looking_time)) + 
  geom_histogram(bins = 90) + 
  xlim(0, 6000)
## Warning: Removed 2 rows containing missing values (geom_bar).

d_no_target %>% 
  ggplot(aes(x = trial_looking_time)) + 
  geom_histogram(bins = 90) + 
  xlim(0, 6000) + 
  facet_wrap(~subject)
## Warning: Removed 84 rows containing missing values (geom_bar).

complexity difference

d %>% 
  ggplot(aes(x = trial_looking_time, 
            fill = trial_stimulus_complexity), 
         ) + 
  geom_density(alpha = 0.5)+ 
  xlim(0, 6000)

d %>% 
  ggplot(aes(x = trial_looking_time, 
            fill = trial_stimulus_complexity), 
         ) + 
  geom_density(alpha = 0.5)+ 
  xlim(0, 6000) + 
  facet_wrap(~subject)

block difference

d %>% 
  ggplot(aes(x = trial_looking_time, 
            fill = block), 
         ) + 
  geom_density(alpha = 0.5)+ 
  xlim(0, 6000)

d %>% 
  ggplot(aes(x = trial_looking_time, 
            fill = block), 
         ) + 
  geom_density(alpha = 0.5)+ 
  xlim(0, 6000) + 
  facet_wrap(~subject)

Complexity differences no target

d_no_target %>% 
  ggplot(aes(x = trial_looking_time, 
            fill = trial_stimulus_complexity), 
         ) + 
  geom_density(alpha = 0.5)+ 
  xlim(0, 6000)

d_no_target %>% 
  ggplot(aes(x = trial_looking_time, 
            fill = trial_stimulus_complexity), 
         ) + 
  geom_density(alpha = 0.5)+ 
  xlim(0, 6000) + 
  facet_wrap(~subject)

Block differences no target

d_no_target %>% 
  ggplot(aes(x = trial_looking_time, 
            fill = block), 
         ) + 
  geom_density(alpha = 0.5)+ 
  xlim(0, 6000)

d_no_target %>% 
  ggplot(aes(x = trial_looking_time, 
            fill = block), 
         ) + 
  geom_density(alpha = 0.5)+ 
  xlim(0, 6000) + 
  facet_wrap(~subject)
## Warning: Groups with fewer than two data points have been dropped.
## Warning in max(ids, na.rm = TRUE): no non-missing arguments to max; returning -
## Inf

Visualizing aggregarded looking time

d_sum_individual <- d %>% 
  group_by(subject, block) %>% 
  summarise(
    mean_lt = mean(trial_looking_time, na.rm = TRUE), 
    sd = sd(trial_looking_time, na.rm = TRUE), 
    n = n(), 
    ci_range_95 = qt(1 - (0.05 / 2), n - 1) * (sd/sqrt(n)), 
    ci_ub = mean_lt + ci_range_95, 
    ci_lb = mean_lt - ci_range_95
  )
## `summarise()` regrouping output by 'subject' (override with `.groups` argument)
d_sum <- d %>% 
  group_by(block) %>% 
  summarise(
    mean_lt = mean(trial_looking_time, na.rm = TRUE), 
    sd = sd(trial_looking_time, na.rm = TRUE), 
    n = n(), 
    ci_range_95 = qt(1 - (0.05 / 2), n - 1) * (sd/sqrt(n)), 
    ci_ub = mean_lt + ci_range_95, 
    ci_lb = mean_lt - ci_range_95
  )
## `summarise()` ungrouping output (override with `.groups` argument)

aggregated

d_sum %>% ggplot(aes(x = block, y = mean_lt)) + 
  geom_pointrange(aes(ymin = ci_lb, ymax = ci_ub)) + 
  theme(axis.text.x = element_text(angle = 45, vjust = 1, hjust=1))

individual

something weird happened

d_sum_individual %>% 
  ggplot(aes(x = block, y = mean_lt)) + 
  geom_pointrange(aes(ymin = ci_lb, ymax = ci_ub)) +
  theme(axis.text.x = element_text(angle = 45, vjust = 1, hjust=1))

d_sum_individual %>% 
  ggplot(aes(x = block, y = mean_lt)) + 
  geom_pointrange(aes(ymin = ci_lb, ymax = ci_ub))  + 
  facet_wrap(~subject) + 
  theme(axis.text.x = element_text(angle = 45, vjust = 1, hjust=1))

weird ones

We can see three very weird ones: SS1604513317537 SS1604515995769 SS1604516882157

weird <- c("SS1604513317537", 
           "SS1604515995769", 
           "SS1604516660396")

d_sum_individual %>% 
  filter(subject %in% weird) %>% 
  ggplot(aes(x = block, y = mean_lt)) + 
  geom_pointrange(aes(ymin = ci_lb, ymax = ci_ub))  + 
  facet_wrap(~subject) + 
  theme(axis.text.x = element_text(angle = 45, vjust = 1, hjust=1))

d %>% 
  filter(subject %in% weird) %>% 
  ggplot(aes(x = trial_looking_time)) + 
  geom_histogram(bins = 90) +
  facet_wrap(~subject) + 
  xlim(0, 6000)+
  theme(axis.text.x = element_text(angle = 45, vjust = 1, hjust=1))
## Warning: Removed 6 rows containing missing values (geom_bar).

Excluding weird ones

d_no_weird <- d %>% 
  filter(!(subject %in% weird)) 

  
d_no_weird_sum <- d_no_weird %>% 
  group_by(block) %>% 
  summarise(
    mean_lt = mean(trial_looking_time, na.rm = TRUE), 
    sd = sd(trial_looking_time, na.rm = TRUE), 
    n = n(), 
    ci_range_95 = qt(1 - (0.05 / 2), n - 1) * (sd/sqrt(n)), 
    ci_ub = mean_lt + ci_range_95, 
    ci_lb = mean_lt - ci_range_95
  )
## `summarise()` ungrouping output (override with `.groups` argument)
d_no_weird_sum %>% ggplot(aes(x = block, y = mean_lt)) + 
  geom_pointrange(aes(ymin = ci_lb, ymax = ci_ub)) + 
  theme(axis.text.x = element_text(angle = 45, vjust = 1, hjust=1))

Demographic

SES

education == “Some high school” ~ 1, education == “High school diploma” ~ 2, education == “Associate Degree/Technical certification” ~ 3, education == “Bachelor’s Degree” ~ 4, education == “Master’s Degree” ~ 5, education == “Doctorate/Professional degree” ~ 6, education == “Other” ~ NA_real_

d %>% 
  distinct(subject, .keep_all = TRUE) %>% 
  ggplot(aes(x = demog_education)) + 
  geom_histogram()
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
## Warning: Removed 4 rows containing non-finite values (stat_bin).

null_m <- lmer(trial_looking_time ~ 1 + (1|subject), 
               data = filter(d, !is.na(demog_education)))

edu_m <- lmer(trial_looking_time ~ demog_education + (1|subject), 
     data = filter(d, !is.na(demog_education)))

summary(edu_m)
## Linear mixed model fit by REML ['lmerMod']
## Formula: trial_looking_time ~ demog_education + (1 | subject)
##    Data: filter(d, !is.na(demog_education))
## 
## REML criterion at convergence: 114420.8
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -1.8912 -0.4691 -0.2950  0.0165  7.5395 
## 
## Random effects:
##  Groups   Name        Variance Std.Dev.
##  subject  (Intercept)  58642   242.2   
##  Residual             331678   575.9   
## Number of obs: 7351, groups:  subject, 38
## 
## Fixed effects:
##                 Estimate Std. Error t value
## (Intercept)       450.69     214.44   2.102
## demog_education    71.82      54.47   1.319
## 
## Correlation of Fixed Effects:
##             (Intr)
## demog_edctn -0.983
anova(null_m, edu_m)
## refitting model(s) with ML (instead of REML)
## Data: filter(d, !is.na(demog_education))
## Models:
## null_m: trial_looking_time ~ 1 + (1 | subject)
## edu_m: trial_looking_time ~ demog_education + (1 | subject)
##        Df    AIC    BIC logLik deviance  Chisq Chi Df Pr(>Chisq)
## null_m  3 114448 114468 -57221   114442                         
## edu_m   4 114448 114475 -57220   114440 1.7923      1     0.1807

Age

d %>% 
  distinct(subject, .keep_all = TRUE) %>% 
  ggplot(aes(x = demog_age)) + 
  geom_histogram()
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
## Warning: Removed 4 rows containing non-finite values (stat_bin).

null_m <- lmer(trial_looking_time ~ 1 + (1|subject), 
               data = filter(d, !is.na(demog_age)))

age_m <- lmer(trial_looking_time ~ demog_age + (1|subject), 
     data = filter(d, !is.na(demog_age)))
## Warning in checkConv(attr(opt, "derivs"), opt$par, ctrl = control$checkConv, :
## Model failed to converge with max|grad| = 0.00430357 (tol = 0.002, component 1)
summary(age_m)
## Linear mixed model fit by REML ['lmerMod']
## Formula: trial_looking_time ~ demog_age + (1 | subject)
##    Data: filter(d, !is.na(demog_age))
## 
## REML criterion at convergence: 114427.7
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -1.8892 -0.4699 -0.2950  0.0175  7.5389 
## 
## Random effects:
##  Groups   Name        Variance Std.Dev.
##  subject  (Intercept)  60451   245.9   
##  Residual             331679   575.9   
## Number of obs: 7351, groups:  subject, 38
## 
## Fixed effects:
##             Estimate Std. Error t value
## (Intercept)  637.516    123.180   5.175
## demog_age      2.458      3.142   0.782
## 
## Correlation of Fixed Effects:
##           (Intr)
## demog_age -0.945
## convergence code: 0
## Model failed to converge with max|grad| = 0.00430357 (tol = 0.002, component 1)
anova(null_m, age_m)
## refitting model(s) with ML (instead of REML)
## Data: filter(d, !is.na(demog_age))
## Models:
## null_m: trial_looking_time ~ 1 + (1 | subject)
## age_m: trial_looking_time ~ demog_age + (1 | subject)
##        Df    AIC    BIC logLik deviance  Chisq Chi Df Pr(>Chisq)
## null_m  3 114448 114468 -57221   114442                         
## age_m   4 114449 114477 -57220   114441 0.6398      1     0.4238

number of background repetition on looking time

full, aggregated, by complexity

full_aggregated <- d %>% 
  mutate(
    number = 1) %>% 
  group_by(
  subject, block, trial_stimulus_type
  ) %>% 
  mutate(num_times_stimulus_seen = cumsum(number))

full_aggregated %>% 
  filter(trial_stimulus_type == "background") %>% 
  ggplot(
    aes(y = log(trial_looking_time), 
        x = num_times_stimulus_seen,
        color = trial_stimulus_complexity)
  ) + 
  geom_point(aes(alpha = 0.2), size = 3, shape = ".") + 
  guides(alpha = FALSE) + 
  labs(color = "Stimulus Complexity") + 
  ylab("Mean Looking Time (log)") + 
  xlab("Number of Stimulus Reptitions") + 
  geom_smooth(method = "lm") + 
  theme(axis.text = element_text(size = 10))
## `geom_smooth()` using formula 'y ~ x'

full, aggregated, by blocks

full_aggregated %>% 
  filter(trial_stimulus_type == "background") %>% 
  ggplot(
    aes(y = log(trial_looking_time), 
        x = num_times_stimulus_seen,
        color = block)
  ) + 
  geom_point(aes(alpha = 0.2), size = 3, shape = ".") + 
  guides(alpha = FALSE) + 
  labs(color = "Stimulus Complexity") + 
  ylab("Mean Looking Time (log)") + 
  xlab("Number of Stimulus Reptitions") + 
  geom_smooth(method = "lm") + 
  theme(axis.text = element_text(size = 10))
## `geom_smooth()` using formula 'y ~ x'

no weird, aggregated, by complexity

excluded_sum <- full_aggregated %>% 
  filter(trial_stimulus_type == "background") %>% 
  filter(!(subject %in% weird))

excluded_sum %>% 
  ggplot(
    aes(y = log(trial_looking_time), 
        x = num_times_stimulus_seen,
        color = trial_stimulus_complexity)
  ) + 
  geom_point(aes(alpha = 0.2), size = 3, shape = ".") + 
  guides(alpha = FALSE) + 
  labs(color = "Stimulus Complexity") + 
  ylab("Mean Looking Time (log)") + 
  xlab("Number of Stimulus Reptitions") + 
  geom_smooth(method = "lm") + 
  theme(axis.text = element_text(size = 10))
## `geom_smooth()` using formula 'y ~ x'

no weird, aggregated, by blocks

excluded_sum %>% 
  ggplot(
    aes(y = log(trial_looking_time), 
        x = num_times_stimulus_seen,
        color = block)
  ) + 
  geom_point(aes(alpha = 0.2), size = 3, shape = ".") + 
  guides(alpha = FALSE) + 
  labs(color = "Stimulus Complexity") + 
  ylab("Mean Looking Time (log)") + 
  xlab("Number of Stimulus Reptitions") + 
  geom_smooth(method = "lm") + 
  theme(axis.text = element_text(size = 10))
## `geom_smooth()` using formula 'y ~ x'

individual, by complexity

full_aggregated %>% 
  filter(trial_stimulus_type == "background") %>% 
  ggplot(
    aes(y = log(trial_looking_time), 
        x = num_times_stimulus_seen,
        color = trial_stimulus_complexity)
  ) + 
  geom_point(aes(alpha = 0.2), size = 3, shape = ".") + 
  guides(alpha = FALSE) + 
  labs(color = "Stimulus Complexity") + 
  ylab("Mean Looking Time (log)") + 
  xlab("Number of Stimulus Reptitions") + 
  geom_smooth(method = "lm") + 
  theme(axis.text = element_text(size = 10)) + 
  facet_wrap(~subject)
## `geom_smooth()` using formula 'y ~ x'

individual, by block

full_aggregated %>% 
  filter(trial_stimulus_type == "background") %>% 
  ggplot(
    aes(y = log(trial_looking_time), 
        x = num_times_stimulus_seen,
        color = block)
  ) + 
  geom_point(aes(alpha = 0.2), size = 3, shape = ".") + 
  guides(alpha = FALSE) + 
  labs(color = "Stimulus Complexity") + 
  ylab("Mean Looking Time (log)") + 
  xlab("Number of Stimulus Reptitions") + 
  geom_smooth(method = "lm") + 
  theme(axis.text = element_text(size = 10)) + 
  facet_wrap(~subject)
## `geom_smooth()` using formula 'y ~ x'

Model

linear

block

null_m <- lmer(trial_looking_time ~ 1 + (1|subject), 
               data = d)

basic_m <- lmer(trial_looking_time ~ block + (1|subject), 
     data = d)

anova(null_m, basic_m)
## refitting model(s) with ML (instead of REML)
## Data: d
## Models:
## null_m: trial_looking_time ~ 1 + (1 | subject)
## basic_m: trial_looking_time ~ block + (1 | subject)
##         Df    AIC    BIC logLik deviance  Chisq Chi Df Pr(>Chisq)   
## null_m   3 121133 121154 -60563   121127                            
## basic_m  6 121126 121167 -60557   121114 13.076      3   0.004474 **
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

stimuli repetition

null_m <- lmer(log(trial_looking_time) ~ 1 + (1|subject), 
               data = excluded_sum)

rep_m <- lmer(log(trial_looking_time) ~ num_times_stimulus_seen + (1|subject), 
     data = excluded_sum)

stimuli repetition and block interaction

interaction_m <- lmer(log(trial_looking_time) ~ num_times_stimulus_seen * block + (1|subject), 
     data = excluded_sum)
anova(rep_m, null_m)
## refitting model(s) with ML (instead of REML)
## Data: excluded_sum
## Models:
## null_m: log(trial_looking_time) ~ 1 + (1 | subject)
## rep_m: log(trial_looking_time) ~ num_times_stimulus_seen + (1 | subject)
##        Df    AIC    BIC  logLik deviance  Chisq Chi Df Pr(>Chisq)    
## null_m  3 6049.9 6069.7 -3022.0   6043.9                             
## rep_m   4 5984.9 6011.3 -2988.5   5976.9 67.016      1  2.693e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
anova(rep_m, interaction_m)
## refitting model(s) with ML (instead of REML)
## Data: excluded_sum
## Models:
## rep_m: log(trial_looking_time) ~ num_times_stimulus_seen + (1 | subject)
## interaction_m: log(trial_looking_time) ~ num_times_stimulus_seen * block + (1 | 
## interaction_m:     subject)
##               Df    AIC    BIC  logLik deviance  Chisq Chi Df Pr(>Chisq)  
## rep_m          4 5984.9 6011.3 -2988.5   5976.9                           
## interaction_m 10 5980.6 6046.5 -2980.3   5960.6 16.322      6    0.01213 *
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
summary(rep_m)
## Linear mixed model fit by REML ['lmerMod']
## Formula: log(trial_looking_time) ~ num_times_stimulus_seen + (1 | subject)
##    Data: excluded_sum
## 
## REML criterion at convergence: 5994.1
## 
## Scaled residuals: 
##     Min      1Q  Median      3Q     Max 
## -3.7972 -0.5324 -0.1798  0.2486  6.3457 
## 
## Random effects:
##  Groups   Name        Variance Std.Dev.
##  subject  (Intercept) 0.1051   0.3241  
##  Residual             0.1726   0.4155  
## Number of obs: 5370, groups:  subject, 39
## 
## Fixed effects:
##                           Estimate Std. Error t value
## (Intercept)              6.3118005  0.0531764 118.695
## num_times_stimulus_seen -0.0046543  0.0005668  -8.211
## 
## Correlation of Fixed Effects:
##             (Intr)
## nm_tms_stm_ -0.189

GAM

only focusing on looking time at the background

library(mgcv)
## Loading required package: nlme
## 
## Attaching package: 'nlme'
## The following object is masked from 'package:lme4':
## 
##     lmList
## The following object is masked from 'package:dplyr':
## 
##     collapse
## This is mgcv 1.8-33. For overview type 'help("mgcv-package")'.
gam_d <- full_aggregated %>% 
  filter(trial_stimulus_type == "background") %>% 
  mutate(trial_stimulus_complexity = as.factor(trial_stimulus_complexity), 
         block = as.factor(block))

number of repetition

gam_m <- gam(trial_looking_time ~ s(num_times_stimulus_seen), 
             data = gam_d, 
             method = "REML")
summary(gam_m)
## 
## Family: gaussian 
## Link function: identity 
## 
## Formula:
## trial_looking_time ~ s(num_times_stimulus_seen)
## 
## Parametric coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  604.987      6.635   91.19   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Approximate significance of smooth terms:
##                              edf Ref.df     F p-value    
## s(num_times_stimulus_seen) 7.879  8.672 21.27  <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## R-sq.(adj) =  0.0324   Deviance explained = 3.38%
## -REML =  42462  Scale est. = 2.4513e+05  n = 5569
plot(gam_m, 
     pages = 1, 
     se = TRUE, 
     shade = TRUE)

gam.check(gam_m)

## 
## Method: REML   Optimizer: outer newton
## full convergence after 6 iterations.
## Gradient range [-0.0001646932,2.934368e-05]
## (score 42462.4 & scale 245132.1).
## Hessian positive definite, eigenvalue range [2.4001,2783.504].
## Model rank =  10 / 10 
## 
## Basis dimension (k) checking results. Low p-value (k-index<1) may
## indicate that k is too low, especially if edf is close to k'.
## 
##                              k'  edf k-index p-value
## s(num_times_stimulus_seen) 9.00 7.88    0.99    0.15
concurvity(gam_m, full = TRUE)
##                  para s(num_times_stimulus_seen)
## worst    5.278689e-27               5.196644e-27
## observed 5.278689e-27               7.734321e-28
## estimate 5.278689e-27               2.303504e-29

number of repetition (with log?)

gam_m <- gam(log(trial_looking_time) ~ s(num_times_stimulus_seen), 
             data = gam_d, 
             method = "REML")
summary(gam_m)
## 
## Family: gaussian 
## Link function: identity 
## 
## Formula:
## log(trial_looking_time) ~ s(num_times_stimulus_seen)
## 
## Parametric coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 6.228538   0.006969   893.8   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Approximate significance of smooth terms:
##                              edf Ref.df     F p-value    
## s(num_times_stimulus_seen) 7.703  8.572 20.92  <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## R-sq.(adj) =  0.0315   Deviance explained = 3.29%
## -REML = 4279.6  Scale est. = 0.27043   n = 5569
plot(gam_m, 
     pages = 1, 
     se = TRUE, 
     shade = TRUE)

gam.check(gam_m)

## 
## Method: REML   Optimizer: outer newton
## full convergence after 6 iterations.
## Gradient range [-0.0003978788,-7.401346e-08]
## (score 4279.63 & scale 0.270434).
## Hessian positive definite, eigenvalue range [2.26755,2783.504].
## Model rank =  10 / 10 
## 
## Basis dimension (k) checking results. Low p-value (k-index<1) may
## indicate that k is too low, especially if edf is close to k'.
## 
##                             k' edf k-index p-value  
## s(num_times_stimulus_seen) 9.0 7.7    0.98   0.045 *
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
concurvity(gam_m, full = TRUE)
##                  para s(num_times_stimulus_seen)
## worst    5.278689e-27               5.196644e-27
## observed 5.278689e-27               5.868023e-28
## estimate 5.278689e-27               2.303504e-29

number of repetition by stimulus complexity

gam_m <- gam(trial_looking_time ~ s(num_times_stimulus_seen, 
                                    by = trial_stimulus_complexity), 
             data = gam_d, 
             method = "REML")

summary(gam_m)
## 
## Family: gaussian 
## Link function: identity 
## 
## Formula:
## trial_looking_time ~ s(num_times_stimulus_seen, by = trial_stimulus_complexity)
## 
## Parametric coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  604.963      6.635   91.18   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Approximate significance of smooth terms:
##                                                               edf Ref.df      F
## s(num_times_stimulus_seen):trial_stimulus_complexitycomplex 6.950  8.031  9.027
## s(num_times_stimulus_seen):trial_stimulus_complexitysimple  7.233  8.253 13.584
##                                                             p-value    
## s(num_times_stimulus_seen):trial_stimulus_complexitycomplex  <2e-16 ***
## s(num_times_stimulus_seen):trial_stimulus_complexitysimple   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## R-sq.(adj) =  0.0325   Deviance explained = 3.49%
## -REML =  42465  Scale est. = 2.4512e+05  n = 5569
plot(gam_m, 
     pages = 1, 
     se = TRUE, 
     shade = TRUE)

gam.check(gam_m)

## 
## Method: REML   Optimizer: outer newton
## full convergence after 7 iterations.
## Gradient range [-0.0001139121,4.48257e-06]
## (score 42465.35 & scale 245121.1).
## Hessian positive definite, eigenvalue range [1.675453,2783.007].
## Model rank =  19 / 19 
## 
## Basis dimension (k) checking results. Low p-value (k-index<1) may
## indicate that k is too low, especially if edf is close to k'.
## 
##                                                               k'  edf k-index
## s(num_times_stimulus_seen):trial_stimulus_complexitycomplex 9.00 6.95    1.02
## s(num_times_stimulus_seen):trial_stimulus_complexitysimple  9.00 7.23    1.02
##                                                             p-value
## s(num_times_stimulus_seen):trial_stimulus_complexitycomplex    0.86
## s(num_times_stimulus_seen):trial_stimulus_complexitysimple     0.86
concurvity(gam_m, full = TRUE)
##                  para
## worst    7.391881e-05
## observed 7.391881e-05
## estimate 7.391881e-05
##          s(num_times_stimulus_seen):trial_stimulus_complexitycomplex
## worst                                                   3.694972e-05
## observed                                                1.413731e-07
## estimate                                                1.013735e-05
##          s(num_times_stimulus_seen):trial_stimulus_complexitysimple
## worst                                                  3.697182e-05
## observed                                               2.975684e-09
## estimate                                               1.001321e-05

number of repetition by stimulus complexity (log?)

gam_m <- gam(log(trial_looking_time) ~ s(num_times_stimulus_seen, 
                                    by = trial_stimulus_complexity), 
             data = gam_d, 
             method = "REML")

summary(gam_m)
## 
## Family: gaussian 
## Link function: identity 
## 
## Formula:
## log(trial_looking_time) ~ s(num_times_stimulus_seen, by = trial_stimulus_complexity)
## 
## Parametric coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 6.228497   0.006967   894.1   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Approximate significance of smooth terms:
##                                                               edf Ref.df      F
## s(num_times_stimulus_seen):trial_stimulus_complexitycomplex 6.517  7.656  7.622
## s(num_times_stimulus_seen):trial_stimulus_complexitysimple  7.146  8.187 15.219
##                                                             p-value    
## s(num_times_stimulus_seen):trial_stimulus_complexitycomplex  <2e-16 ***
## s(num_times_stimulus_seen):trial_stimulus_complexitysimple   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## R-sq.(adj) =  0.0322   Deviance explained = 3.45%
## -REML = 4287.1  Scale est. = 0.27026   n = 5569
plot(gam_m, 
     pages = 1, 
     se = TRUE, 
     shade = TRUE)

gam.check(gam_m)

## 
## Method: REML   Optimizer: outer newton
## full convergence after 6 iterations.
## Gradient range [-0.0004995101,2.60902e-05]
## (score 4287.147 & scale 0.2702625).
## Hessian positive definite, eigenvalue range [1.622075,2783.007].
## Model rank =  19 / 19 
## 
## Basis dimension (k) checking results. Low p-value (k-index<1) may
## indicate that k is too low, especially if edf is close to k'.
## 
##                                                               k'  edf k-index
## s(num_times_stimulus_seen):trial_stimulus_complexitycomplex 9.00 6.52    1.02
## s(num_times_stimulus_seen):trial_stimulus_complexitysimple  9.00 7.15    1.02
##                                                             p-value
## s(num_times_stimulus_seen):trial_stimulus_complexitycomplex    0.90
## s(num_times_stimulus_seen):trial_stimulus_complexitysimple     0.88
concurvity(gam_m, full = TRUE)
##                  para
## worst    7.391881e-05
## observed 7.391881e-05
## estimate 7.391881e-05
##          s(num_times_stimulus_seen):trial_stimulus_complexitycomplex
## worst                                                   3.694972e-05
## observed                                                1.021689e-07
## estimate                                                1.013735e-05
##          s(num_times_stimulus_seen):trial_stimulus_complexitysimple
## worst                                                  3.697182e-05
## observed                                               1.087388e-07
## estimate                                               1.001321e-05

number of repetition by block

gam_m <- gam(trial_looking_time ~ s(num_times_stimulus_seen, by = block), 
             data = gam_d, 
             method = "REML")

plot(gam_m, all.terms = TRUE, pages = 1,
     shade = TRUE,
     rug = TRUE, 
     se = TRUE)

summary(gam_m)
## 
## Family: gaussian 
## Link function: identity 
## 
## Formula:
## trial_looking_time ~ s(num_times_stimulus_seen, by = block)
## 
## Parametric coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  604.966      6.632   91.22   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Approximate significance of smooth terms:
##                                                         edf Ref.df      F
## s(num_times_stimulus_seen):blockall_complex           5.876  7.035  4.642
## s(num_times_stimulus_seen):blockall_simple            7.234  8.255 10.986
## s(num_times_stimulus_seen):blockmixed_complex_deviant 4.119  5.080  6.754
## s(num_times_stimulus_seen):blockmixed_simple_deviant  5.588  6.735  4.952
##                                                        p-value    
## s(num_times_stimulus_seen):blockall_complex           3.15e-05 ***
## s(num_times_stimulus_seen):blockall_simple             < 2e-16 ***
## s(num_times_stimulus_seen):blockmixed_complex_deviant 2.17e-06 ***
## s(num_times_stimulus_seen):blockmixed_simple_deviant  2.04e-05 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## R-sq.(adj) =  0.0333   Deviance explained = 3.72%
## -REML =  42462  Scale est. = 2.4492e+05  n = 5569
gam.check(gam_m)

## 
## Method: REML   Optimizer: outer newton
## full convergence after 8 iterations.
## Gradient range [-0.0004214249,4.232332e-05]
## (score 42461.66 & scale 244917.5).
## Hessian positive definite, eigenvalue range [0.9508617,2782.009].
## Model rank =  37 / 37 
## 
## Basis dimension (k) checking results. Low p-value (k-index<1) may
## indicate that k is too low, especially if edf is close to k'.
## 
##                                                         k'  edf k-index p-value
## s(num_times_stimulus_seen):blockall_complex           9.00 5.88       1    0.43
## s(num_times_stimulus_seen):blockall_simple            9.00 7.23       1    0.44
## s(num_times_stimulus_seen):blockmixed_complex_deviant 9.00 4.12       1    0.36
## s(num_times_stimulus_seen):blockmixed_simple_deviant  9.00 5.59       1    0.33
concurvity(gam_m, full = TRUE)
##                  para s(num_times_stimulus_seen):blockall_complex
## worst    0.0001035554                                2.539898e-05
## observed 0.0001035554                                3.225119e-06
## estimate 0.0001035554                                1.611805e-05
##          s(num_times_stimulus_seen):blockall_simple
## worst                                  2.329026e-05
## observed                               5.265772e-09
## estimate                               4.459062e-06
##          s(num_times_stimulus_seen):blockmixed_complex_deviant
## worst                                             1.697086e-05
## observed                                          1.606767e-07
## estimate                                          5.689955e-06
##          s(num_times_stimulus_seen):blockmixed_simple_deviant
## worst                                            3.790317e-05
## observed                                         6.899570e-06
## estimate                                         2.979582e-06

number of repetition by block (log?)

gam_m <- gam(log(trial_looking_time) ~ s(num_times_stimulus_seen, by = block), 
             data = gam_d, 
             method = "REML")

plot(gam_m, all.terms = TRUE, pages = 1,
     shade = TRUE,
     rug = TRUE, 
     se = TRUE)

summary(gam_m)
## 
## Family: gaussian 
## Link function: identity 
## 
## Formula:
## log(trial_looking_time) ~ s(num_times_stimulus_seen, by = block)
## 
## Parametric coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  6.22852    0.00697   893.6   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Approximate significance of smooth terms:
##                                                         edf Ref.df      F
## s(num_times_stimulus_seen):blockall_complex           4.655  5.705  4.434
## s(num_times_stimulus_seen):blockall_simple            6.808  7.913 11.433
## s(num_times_stimulus_seen):blockmixed_complex_deviant 4.471  5.493  7.375
## s(num_times_stimulus_seen):blockmixed_simple_deviant  5.068  6.171  3.938
##                                                        p-value    
## s(num_times_stimulus_seen):blockall_complex           0.000393 ***
## s(num_times_stimulus_seen):blockall_simple             < 2e-16 ***
## s(num_times_stimulus_seen):blockmixed_complex_deviant  1.3e-06 ***
## s(num_times_stimulus_seen):blockmixed_simple_deviant  0.000487 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## R-sq.(adj) =  0.0312   Deviance explained = 3.49%
## -REML = 4299.7  Scale est. = 0.27053   n = 5569
gam.check(gam_m)

## 
## Method: REML   Optimizer: outer newton
## full convergence after 7 iterations.
## Gradient range [-0.0005736542,3.267841e-05]
## (score 4299.665 & scale 0.2705294).
## Hessian positive definite, eigenvalue range [0.486995,2782.007].
## Model rank =  37 / 37 
## 
## Basis dimension (k) checking results. Low p-value (k-index<1) may
## indicate that k is too low, especially if edf is close to k'.
## 
##                                                         k'  edf k-index p-value
## s(num_times_stimulus_seen):blockall_complex           9.00 4.66    0.99    0.34
## s(num_times_stimulus_seen):blockall_simple            9.00 6.81    0.99    0.27
## s(num_times_stimulus_seen):blockmixed_complex_deviant 9.00 4.47    0.99    0.26
## s(num_times_stimulus_seen):blockmixed_simple_deviant  9.00 5.07    0.99    0.31
concurvity(gam_m, full = TRUE)
##                  para s(num_times_stimulus_seen):blockall_complex
## worst    0.0001035554                                2.539898e-05
## observed 0.0001035554                                6.175376e-06
## estimate 0.0001035554                                1.611805e-05
##          s(num_times_stimulus_seen):blockall_simple
## worst                                  2.329026e-05
## observed                               2.196112e-08
## estimate                               4.459062e-06
##          s(num_times_stimulus_seen):blockmixed_complex_deviant
## worst                                             1.697086e-05
## observed                                          4.328620e-10
## estimate                                          5.689955e-06
##          s(num_times_stimulus_seen):blockmixed_simple_deviant
## worst                                            3.790317e-05
## observed                                         7.738447e-06
## estimate                                         2.979582e-06

target reaction time predictive of anything?

import the complexity norm